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ABSTRACT 

Two-equation turbulence models are analyzed from the 
perspective of spectral closure theories. Kolmogorov the- 
ory provides useful information for models, but it is lim- 
ited to ‘equilibrium’ conditions in which the energy spectrum 
has relaxed to a steady state consistent with the forcing at 
large scales; it does not describe transient evolution between 
such states. Transient evolution is necessarily through ‘non- 
equilibrium’ states, which can only be found from a theory of 
turbulence evolution, such as one provided by a spectral clo- 
sure. When the departure from equilibrium is small, perturba- 
tion theory can be used to approximate the evolution by a two- 
equation model. The perturbation theory also gives explicit 
conditions under which this model can be valid, and when it 
will fail. Implications of the non-equilibrium corrections for 
the classic Tennekes-Lumley balance in the dissipation rate 
equation are drawn: it is possible to establish both the can- 
cellation of the leading order Re 1 / 2 divergent contributions to 
vortex stretching and enstrophy destruction, and the existence 
of a nonzero difference which is finite in the limit of infinite 
Reynolds number. 


INTRODUCTION 

The two-equation model is solved as part of many 
Reynolds averaged models, even in models that include 
Reynolds stresses and other more refined descriptions of the 
large-scale statistics. By itself, a two-equation model treats 
turbulent fluctuations as homogeneous isotropic turbulence. 
This simplification can be motivated by appealing to Kol- 
mogorov theory, which supports the translation and rotation 
invariance of the statistics of fluctuations at small scales, re- 
gardless of the inhomogeneity and anisotropy of the large 
scales (Frisch, 1995). This paper analyzes the two-equation 
model as a description of the evolution of statistically non- 
stationary homogeneous isotropic turbulence. It poses the 
question of when two-equation modeling is possible, and 
when is it impossible. The answer is closely linked to the 
classic Tennekes and Lumley (1972) balance between vortex 
stretching and enstrophy destruction. 


NON-STATIONARY TURBULENCE 

Kolmogorov theory implies that statistically stationary 
homogeneous isotropic turbulence is characterized by two pa- 
rameters, equivalent to a length scale £ and the dissipation rate 
e (this choice of variables will be justified later). This idea is 
formulated in terms of the energy spectrum E ( k) as 

E(k)=E{k;£,e) (1) 

Such a state is sometimes called an ‘equilibrium.’ but it is 
more correct to call it simply a steady state, since thermody- 
namically, it is far from equilibrium. Eq. (1) includes a non- 
trivial idea from Kolmogorov’s theory, that because the dissi- 
pation rate equals the energy flux through the inertial range, it 
links the large and small scales; this justifies the characteriza- 
tion of the entire energy spectrum by a flux, e, and a largest 
dynamically relevant scale, £. 

But these considerations are purely ‘static;’ modeling 
attempts to treat temporal evolution under unsteady condi- 
tions. For the purposes of modeling, a non-stationary, or time- 
dependent energy spectrum cannot vary arbitrarily in time; for 
a finite dimensional model to be possible at all, it can only 
depend on time through finitely many parameters, which can 
reasonably be taken to be l and £ themselves. Then Eq. (1) is 
replaced by the time-dependent ansatz 

E(K,t)=E(K-J(t),e(t)) (2) 

This might be called a (temporally) local ‘equilibrium’ or 
steady state. This special solution structure, in which the time 
dependence occurs only through the time dependence of two 
parameters, can be compared to the normal solution of kinetic 
theory (Chapman and Cowling, 1970). 

But how can evolution equations for £ and e be found? 
Kolmogorov theory is silent on this point. It states only that 
that starting from any state, under steady forcing, turbulence 
will evolve to an equilibrium steady state with l and e deter- 
mined by the forcing. Thus, Kolmogorov theory merely as- 
serts the existence and analytical form of an attractor for the 
statistics of small scales of motion; it provides no information 
about the dynamics of reaching that state. 
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The usual approach to modeling postulates two equa- 
tions: one is an exact energy balance, but the second is basi- 
cally a fabrication. Briefly, it states that in a self-similar flow, 
since all quantities with the same units are proportional, the 
frequencies e/e, e/k, and P/k ( k is the turbulent kinetic en- 
ergy and P is the rate of energy input) are constant multiples 
of each other, so that 



d K 4VK^(K,t) 


and enstrophy destruction (or palinstrophy ) G is defined by 


G(t)= [ dK4v 2 K 4 E(K,t). (8) 

Jo 


e P e 

e= C *'k- C *k (3) 

The constants are determined by comparison with numerical 
or experimental data in self-similar flows, usually self-similar 
decaying turbulence and homogeneous shear flow, in which, 
to repeat, a model of this form is necessarily correct because 
of self-similarity (Clark. 1999). The extension of this equation 
to general flows assumes a ‘continuum mechanics’ paradigm, 
in which turbulence is viewed as a ‘fluid’ with the same 'prop- 
erties’ in every flow. However, this interesting paradigm is 
introduced without fundamental justification. 

It is unnecessary to insist that this approach is scientifi- 
cally unsatisfactory, regardless of any practical value it might 
be found to have, and indeed, this fact is well recognized in the 
modeling literature. Improvement requires equations of mo- 
tion for E{K,t). It follows from the Navier-Stokes equations 
and the kinematic assumptions of homogeneity and isotropy 
that 

d , 

E(K,t) ~n(K,t)~ -=-&(K,t) — 2vK-E(K,t) (4) 

u K 


The integration by parts in Eq. (7) is justified by the plausi- 
ble assumption of exponential decay of turbulence statistics 
when K — > oo, so that \im K ^oo2v K 2 ^(K,t) = 0. Since pro- 
duction is assumed concentrated at large scales, the integral 

/ dK 2 vk 2 H (kj) is of order Re~ ] and is ignored in Eq. 
JO 

(6) in the high Reynolds number limit. These results can also 
be obtained from the Karman-Howarth equation (Batchelor, 
1953). 


EQUILIBRIUM EVOLUTION 

Important properties of an equilibrium state in which J? 
is constant and E ~ K"~ 5 /3 in the inertial range follow because 
the integrals defining 5 and G, Eqs. (7) and (8), are dominated 
by contributions near the Kolmogorov scale Kj <x (e/v 3 ) 1 / 4 
(Smith and Reynolds, 1991). Then 

S ss vicj « v^ 1 / 2 ~ Re 1 / 2 (9) 

and 


where n is a known, possibly time-dependent forcing (or 
‘production’) spectrum, and where J? is the energy flux, a 
third order moment of the velocity field (Batchelor, 1953). 

An energy balance is obtained by integrating Eq. (4) over 
all wavenumbers K, 

k(t) = P(t)-£(t) (5) 

where total kinetic energy k(t) = / dK E(K,t), total produc- 

Jo 

tion rate P(t) = / dK H(K,t), and dissipation rate e(t) = 

Jo 

/ dK 2vK 2 E(K,t). The flux term integrates to zero because 

Jo 

it is a gradient; recall that the gradient form reflects the con- 
servative nature of energy transfer by nonlinear interaction. 
While this equation is attractive because it is exact, it must 
be emphasized that the nonlinearity embodied in the flux term 
J? , the heart of the turbulence problem, plays no role in this 
equation. 

Similarly, a dissipation rate equation is obtained by mul- 
tiplying Eq. (4) by 2 vie 2 and integrating, 

£ = S(t) - G(t) (6) 


G ss v 2 iej 0/3 oe v~ 1//2 ~ Re 1//2 , (10) 

the famous Re 1 / 2 divergences identified by Tennekes and 
Lumley (1972). Since in an equilibrium state, e = S — G = 0, 
it follows that S = G also in a local equilibrium defined by 
Eq. (2). Therefore, evolution through exactly local equilib- 
rium states is inconsistent with evolution of £: any attempt to 
derive an e equation based on the assumption of a local equi- 
librium must be inconclusive. 

Further consequences of equilibrium evolution follow 
front the Kolmogorov similarity form of the energy spectrum 
for steady state turbulence 

E{K-e,£) = £ 2 ^£ 5 ^(j>(eK) (ii) 

The function <j> is ‘universal’ insofar as it does not depend 
on £ or e, however, the particular analytical form of 0 only 
applies in one class of flows: forced steady state homogeneous 
isotropic turbulence. From 

k = e 2 / 3 f/ 3 f dx(j)(x) ( 12 ) 

Jo 


where vortex stretching S is defined by 


the Kolmogorov relation 


S(t) 



dK 2VK 2 —^ r (K,t) 
dK 


£ = C e 


k 2 / 2 

£ 


(13) 
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with 



follows. The constant C e is ‘universal’ in the restricted sense 
noted above: whereas Eq. (13) can apply more generally, for 
example to decaying turbulence, the dimensionless spectrum 
0 can be flow-dependent, so that the constant C £ can take dif- 
ferent values in different flows (Bos et al., 2007). 

Similar considerations apply to the spectral eddy viscos- 
ity, which will be discussed further in the next section, 


v, 



Inserting the self-similar form Eq. (11), 


v, = 


J dx y 0(x)/.r 3 k 2 



(15) 


(16) 


The spectral eddy viscosity is therefore a constant multiple 
of the two-equation model viscosity Ar 2 /e; the proportionality 
constant has a theoretical value that does not depend on £ and 
£, but that could depend on the flow class through the function 
0 (Besnard et al, 1996). 

Because they refer to equilibrium conditions alone, these 
results do not depend on any particular expression for the flux 
P in Eq. (4). This fact can account for their ubiquity in tur- 
bulence modeling. 


NON-EQUILIBRIUM EVOLUTION 

While considerations based on Kolmogorov theory are 
very useful, they only carry us so far. Woodruff and Ru- 
binstein (2006) advocated a new viewpoint on two-equation 
modeling, that the two-equation model should determine an 
approximate solution of the spectral evolution equation Eq. 
(4) under some closure assumption for the energy flux P . For 
analytical simplicity. Woodruff and Rubinstein (2006) used 
the classical Heisenberg model (Batchelor, 1953), in which 
the energy flux P is the functional of the energy spectrum. 


model to more complete closure theories based on triad inter- 
actions (Kraichnan, 1987; Rubinstein and Clark, 2004). The 
Heisenberg model treats the energy flux at wavenumber k: as 
the product of a scale-dependent eddy viscosity at wavenum- 
ber K, the second integral on the right side of Eq. (17), and 
the square of the effective strain at wavenumber K , given by 
the first integral on the right side of Eq. (17). 

The time-independent form of Eq. (4), 

d , 

0 = P(k)- -^P(k) -2vk~E(k) (19) 


has a steady solution consistent with the Kolmogorov K"' 5 / 3 
inertial range; then the constant C /, can be chosen to match 
the familiar Kolmogorov constant, treated as an empirical in- 
put. The existence of a Kolmogorov steady state solution es- 
tablishes the consistency of the Heisenberg model with Kol- 
mogorov theory, but the crucial fact is that the Heisenberg 
model is more general, because it states equations for the evo- 
lution of the energy spectrum without imposing Kolmogorov 
scaling in advance. From this viewpoint, spectral closures 
like the Heisenberg model stand in the same relation to Kol- 
mogorov theory as the Boltzmann equation to equilibrium 
thermodynamic s . 

The forcing spectrum n()C,f) is assumed to be com- 
pletely specified by a single length-scale L(t) such that n( K, t ) 
is nonzero only near the wavenumber k/L, and by the total 

energy input rate P(t) = f dK n(ic,f). If P and L do not 


depend on time, then the corresponding solution of Eq. (19) 
is 


E eq (K) = E(K-J,e) t = L e = P (20) 

This connection between the energy spectrum and the forcing 
spectrum justifies the somewhat unconventional choice of l 
and e as basic variables. It is obvious that if L = L(t) and 
P = P(t) depend on time, then the local equilibrium 


E lo c(K,t)=E(K-J{t),£(t)) (21) 


with 


m=L(t) e(t) = P{t) (22) 


P[E(k)] = C h [ K dn p 2 E(p) r dp E(p)G(p) (17) 
JO Jk 

with the ‘algebraic’ time- scale 

e[p) =7wm 

The Heisenberg model is written in the form of two equations, 
Eqs. (17) and ( 1 8) in order to emphasize the connection of this 


satisfies 

d , 

0 = P(K,t)- ^-P(K,t)-2VK 2 E(K,t) (23) 

dK 

because the flux function P depends on time only through 
E. But the complete spectral evolution model Eq. (4) is not 
satisfied by E\ 0C (K,t) because of the time derivatives on the 
left. Thus, the spectral evolution equation cannot be satis- 
fied by evolution through local equilibrium states. If the time 
derivatives are simply ignored, then there is no need for a two- 
equation model, because the equalities in Eq. (22) state that 
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the spectrum adjusts instantaneously to the forcing. But if 
the time derivatives are small (this assumption can be made 
precise by perturbation theory: compare (Woodruff and Ru- 
binstein, 2006)), then it is reasonable to construct an approxi- 
mate solution of Eq. (4) by treating the local equilibrium as a 
leading order solution, say 

E 0 (K-J(t),e(t))=E loc (K,t) (24) 

where analytical form in Eq. (21) is assumed, but Eq. (22) is 
not imposed. A perturbative solution of Eq. (4) is then sought 
in the form 

E{K,t) = E 0 {lcJ{t),E{t))+Ei(K,t) (25) 

where £j <C Eq (again, this statement can be made precise 
in perturbation theory). This approach to solving the spectral 
evolution equation is motivated by the Chapman-Enskog and 
Hilbert expansions of kinetic theory (Chapman and Cowling, 
1970). 

The correction E\ satisfies the linear integral equation 

~S[Ei(k)] =-j t E 0 (K-,e(t),e(t)) (26) 

where is the flux linearized about the zero order solution; 
for the Heisenberg model Eq. (17), 

2zf[£i] = C,, / dp p 2 E\{p) f dp 
JO Jk 

+C4 /"V p 2 E 0 (p) [ °° dp 4= (2V) 

2 70 Jk y/p3£ 0 (p) 

Thus. E i cancels the error in Eq. (4) due to the time depen- 
dence of Eq. 

The crux of the analysis of Woodruff and Rubinstein 
(2006) is that the possibility of solving Eq. (26) requires that 
certain compatibility conditions be satisfied, and these condi- 
tions provide the two-equation model for the unknown func- 
tions f(t) and e(t). 

Assuming that the solution exists, it is clear from Eqs. 
(26) and (27) that the correction £j has the structure 

E.IM-AWgj+SMlj <28) 

and Eq. (26) shows that in the inertial range where Eq 

JC- 5 / 3 , 

El cx Kf~ 7/3 (29) 

The existence of such a ‘non-equilibrium' correction to the 
enegy spectrum was first proposed by Yoshizawa (1994). Al- 
though it is very difficult to extract this contribution even in 


closure calculations, by using a novel method to analyze the 
data, Horiuti and Ozawa (2011) have recently given conclu- 
sive evidence for its existence in direct numerical simulations. 

This analysis demonstrates that in a regime described by 
a two-equation model, spectral evolution cannot be through 
exactly equilibrium states, and that non-equilibrium correc- 
tions for the time dependence necessarily arise. These correc- 
tions have a crucial role in answering the famous and funda- 
mental question posed by Tennekes and Lumley (1972): given 
that 5 and G in Eq. (6) are shown in Eqs. (7) and (8) to diverge 
as v- 1 / 2 ~ Re 1 ' 2 , how is it possible for e = 5— G to be finite 
and nonzero as Re —> °o? 

The statement that the difference between vortex stretch- 
ing and enstrophy destruction is independent of Reynolds 
number expresses the crucial property that the evolution of the 
small scales depends entirely on the large scales; or, equiva- 
lently, that e does indeed link the large and small scales, even 
in a time-dependent turbulent flow . This conclusion is abso- 
lutely critical to the possibility of modeling: it bears repeating 
that this conclusion does not follow from Kolmogorov theory. 

S and G can be written in terms of contributions from Eq 
and Ei as 

S = S Q + Si G = Gq + Gi (30) 

but since in local equilibrium, Sq = Go, 

e=S-G = S 1 -G 1 (31) 

To evaluate Si, note that the contribution to the flux from the 
correction £i is 

&l =2zf[£ 1 ] ~ ic~ 2/3 (32) 

in the inertial range, where Ji?[E i] is defined in Eq. (27); con- 
sequently, 

S! ~ J dK vk^i ~ vxf ~ £ ~ Re 0 (33) 

Similarly, 

Gi ~ J dK v 2 K* 4 £j ~ v 2 )^ 3 ~ £ 2 ~ Ee° (34) 

Thus, Si — Gi can be nonzero and finite in the high Reynolds 
number limit. This argument establishes the Tennekes- 
Lumley balance in a perturbation theory of slow spectral evo- 
lution. 

The Tennekes-Lumley balance has a crucial role in the 
formulation of a dissipation rate equation: it states that the 
quantity that must be modeled is the unknown 0(Re a ) differ- 
ence between two correlations of order 0(Re 1 / 2 ); the diffi- 
culty is that this difference does not correspond to any defi- 
nite correlation of flow quantities. Consequently, the program 
of ‘correlation by correlation' modeling is not possible. The 
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result is the usual heuristic approach to the dissipation rate 
equation. The present analysis provides a possible explana- 
tion for this impasse: the dissipation rate equation arises as a 
first order effect in perturbation theory. Since it is therefore 
tightly linked to non-equilibrium effects, it is not accessible to 
derivations based on equilibrium considerations alone. 

The time dependence also induces corrections to the 
equilibrium eddy viscosity and length scale relations in Eqs. 
(13) and (15). The correction to eddy viscosity has the form 


1 /•“ Ei k) 

Vrl = x / dK 

2 do ypk : 3 E q (k) 

i M r dK B w 

2e{t)J 0 y/K^j 


r,- 

2 e(t) Jo s/k^M*) 
(35) 


where Eq. (28) has been substituted for E \ . Non-equilibrium 
corrections to the eddy viscosity have been introduced by 
Yoshizawa and Nisizima (1993). However, their goal was to 
find a correction to the eddy viscosity in steady state turbu- 
lence that makes it appropriate for homogeneous shear flow 
as well. The present perspective is that the eddy viscosity co- 
efficient in homogeneous shear flow will be different from that 
in steady state turbulence because the function <j> in Eq. (16) 
will not be the same. The non-equilibrium corrections in Eq. 
(35) are small corrections to the equilibrium eddy viscosity 
in weakly time-dependent homogeneous isotropic turbulence; 
they do not apply to any other class of turbulent flows. 

The analogous corrections to Eq. (13) are non- 
equilibrium corrections to the Kolmogorov law e ~ A 3 / 2 /£. 
Attempts to derive a dissipation rate equation by manipula- 
tion of this equilibrium law are futile because the dissipation 
rate equation depends on non-equilibrium effects. 


LIMITATIONS OF TWO-EQUATION MODELS 

This analysis demonstrates that two-equation modeling 
requires an explicit treatment of the departure of a turbu- 
lent flow from equilibrium conditions. But it also restricts 
the validity of models to weak departures from equilibrium, 
in which the ratios e/e and £/£, which quantify the depar- 
ture from equilibrium, are ‘small' in a sense made precise 
in the theory. An extreme case of violation of this restric- 
tion occurs if either e/e or £/£ is large enough to cause 
Eq + E\ to become negative. This possibility can be compared 
to ‘continuum breakdown' (Bird, 1994) in kinetic theory, in 
which the Navier-Stokes equations become invalid because 
the Chapman-Enskog distribution is non-realizable. The pos- 
sibility of formulating this limitation depends on a key feature 
of the present analysis, that it constructs an approximate en- 
ergy spectrum, as well as equations of evolution for the turbu- 
lence descriptors e and £. 

It should be stressed that the explicit formulation of con- 
ditions under which the two-equation model can be valid is a 
particular strength of the present perturbative approach; noth- 
ing comparable is possible if the model is simply postulated 
by the usual heuristic arguments. 

Weak departure from equilibrium imposes relations be- 
tween large and small scales: they are linked to leading or- 
der in perturbation theory through the dissipation rate e as if 


the flow were steady, and the evolution is therefore such that 
the small scale evolution is basically determined completely 
by the large scales. This relation is sometimes expressed by 
saying that the small scales are ‘slaved’ to the large scales; 
however, from the perspective of closure theories, this rela- 
tion between the small and large scales is very special, not a 
general feature of all turbulent flows. 

Of course, the two-equation model can always be solved, 
even when it is not valid. In such cases, the relation between 
small and large scales that it imposes can make the evolution 
unrealistic. An example occurs in the problem of homoge- 
neous isotropic turbulence with linearly increasing production 
rate (Rubinstein et al., 2004), in which turbulence evolves 
from a Kolmogorov steady state to a self-similarly growing 
state. A two-equation model consistent with both the initial 
steady state and with the state of self-similar growth misrep- 
resents the evolution of the dissipation rate: whereas closure 
and direct numerical simulation data suggest a delay in the 
growth of e followed by a very rapid rise, two-equation mod- 
els predict a gradual growth of e throughout, as if the true 
dynamics were subjected to a low-pass filter. This ‘filtering’ 
effect arises because the two-equation model links the evolu- 
tion of the dissipation rate to large scales at all times. 

During the initial transient, the large scales are evolving 
more rapidly than the small scales; but when the small scales 
finally start to respond to the large scales, they do so on their 
own fast time scales, not on the time scales of the large scale 
motion. In both cases, the link between large and small scales 
is broken, and the the energy spectrum simply does not admit 
the two-parameter ansatz of Eq. (2). A more familiar example 
of breaking the link between large and small scales is the rise 
in dissipation rate during the initial transient of decaying tur- 
bulence, before the Kolmogorov spectrum is formed. In this 
problem, there is no ‘production’ mechanism at large scales to 
cause the dissipation rate to increase; the small scales evolve 
independently of the large scales, and consequently, the two- 
equation model cannot predict the evolution. Compare also 
the analysis of front propagation in transient turbulent evolu- 
tion by Connaughton and Nazarenko (2004). 

The breakdown of the two-equation model in these prob- 
lems does not imply that a statistical approach is impossible: 
as noted earlier, spectral closure theories do not make any spe- 
cial assumptions about Kolmogorov scaling, and are therefore 
candidates to analyze flows in which a Kolmogorov spectrum 
does not exist. 


EXTENSIONS OF THE THEORY 

The term ‘equilibrium’ is used in turbulence modeling 
very generally, to describe self-similar states, even when they 
are time-dependent, as in decaying turbulence and in homoge- 
neous shear. Another important feature of the present analysis 
is that perturbations of each such ‘equilibrium’ state requires a 
separate analysis: there is no assertion that a model developed 
for departures from steady-state turbulence will correctly de- 
scribe turbulent decay. A possible topic for future research is 
the development along the present lines of a model for pertur- 
bations of self-similar turbulent boundary layers. 
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CONCLUSIONS 

Time-dependent evolution through equilibrium states is 
not possible; any description of statistically unsteady dynam- 
ics requires some theory to quantify the departures from equi- 
librium. When the time dependence is sufficiently slow, per- 
turbation theory can be used to derive a two-equation model 
which provides an approximate description of the evolution. 
The perturbation theory also indicates the limitations of the 
two-equation model. Completely general transient evolution 
can only be predicted by solving a complete spectral evolution 
equation. The use in this analysis of the Heisenberg model 
was entirely for analytical convenience; a complete perturba- 
tive development of turbulence models must be based on more 
accurate closures that include the dynamics of triad interac- 
tions in turbulence (Kraichnan, 1959). 
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